Dipolar Bose gases: Many-body versus mean-field description 
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We characterize zero-temperature dipolar Bose gases under external spherical confinement as a 
function of the dipole strength using the essentially exact many-body diffusion Monte Carlo (DMC) 
technique. We show that the DMC energies are reproduced accurately within a mean-field framework 
if the variation of the s-wave scattering length with the dipole strength is accounted for properly. 
Our calculations suggest stability diagrams and collapse mechanisms of dipolar Bose gases that 
differ significantly from those previously proposed in the literature. 

PACS numbers: 



The recently achieved Bose-Einstein condensation of 
atomic chromium 1, 2] has added two new twists to the 
study of ultracold matter. First, Cr condensates realize 
the first spin-three spinor condensate [j,|a]. Second, they 
exhibit, due to Cr's comparatively large magnetic dipole 
moment, observable anisotropic long-range interactions 
|3|. These long-range interactions allow the relative ori- 
entation between well separated atoms or molecules to 
be controlled, either by tuning external fields, or else by 
adjusting trap anisotropy. An extensive theoretical lit- 
erature has predicted novel properties for these gases. 
For example, Roton-like features have been predicted for 
trapped gases JJ_, J^J, along with unique phases such 
as checkerboard, Mott insulator, and supersolid phases 

Ellllii. 

Rapid experimental progress in cooling and trapping 
suggests that condensation of ground state molecules 



with large permanent dipole moments, such as OH 
RbCs jSj J KRb ,2i] , NH jl^ , may be achieved soon. These 
species would represent truly strongly interacting dipoles, 
with interparticle interaction strengths up to ~ 10^ times 
stronger than in chromium. Indeed, the dipolar interac- 
tions could become the dominant energy scale in such 
systems, driving transitions to strongly correlated states 
of these gases. 

Thus far, dipolar Bose gases at zero temperature have 
been described by the mean- field Gross-Pitaevskii (GP) 
equation. In particular, stability diagrams and excita- 
tion spectra have been derived within this formalism 
lUIIillllllliSIIllilil. Somewhat surprisingly, the 
validity of the GP equation for dipolar gases with strong, 
anisotropic long-range interactions has not been assessed 
in detail to date. Per se it is not clear that a Hartree 
wave function, as used in the GP framework, can prop- 
erly describe systems interacting through potentials that 
fall off as zLl/r^ at large interparticle distances. Neutral 
atom-atom interactions, e.g., fall off as — 1/r^ and mean- 
field treatments are shown to predict the properties of di- 
lute atomic Bose gases with high accuracy. This requires, 
however, replacing the true interaction by an appropri- 



ate Fermi pseudopotential. For electronic systems with a 
repulsive 1/r interactions, a Hartree-Fock formalism is a 
suitable starting point for computing the electronic struc- 
ture of atoms and molecules. An accurate determination 
of observables, however, often requires correlation effects 
beyond those described by a Hartree-Fock wave function. 

To address this issue for dipolar interactions, this Let- 
ter reports essentially exact many-body diffusion Monte 
Carlo (DMC) calculations for dipolar Bose gases inter- 
acting through realistic two-body model potentials. Our 
DMC results indicate that the GP equation is adequate 
to describe the gas, provided that the pseudo-potential is 
parameterized in terms of a "dipole-normalized" s-wave 
scattering length a{d). Using this renormalized a{d) in- 
stead of the "bare" s-wave scattering length suggests dis- 
tinctly different collapse behaviors and stability diagrams 
than proposed in the literature and has important impli- 
cations for the experimental realization of dipolar Bose 
gases. 

Consider the Hamiltonian H for N interacting bosonic 
dipoles with mass m, assumed to be polarized along the 
z-axis, under external harmonic confinement. 
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where uj denotes the trapping frequency, fj the position 
vector with respect to the trap center of the jth dipole, 
and fjk the distance vector, fjk = fj — fk- We model the 
boson-boson potential V{f) by a short-range hardcore 
with cutoff radius b and a long-range tail with dipole 
moment d. 



V{f) = 
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where 9 denotes the angle between the vector f and the 
laboratory z-axis. The length D* = md^/H^, at which 
the characteristic two-body potential and kinetic ener- 
gies coincide, is used in the following to characterize the 
anisotropic long-range interaction. 



The Hamiltonian H in Eq. (^ applies, in cgs units, to 
bosons with either magnetic or electric dipole moments. 
Importantly, the induced dipole moments that drive the 
interaction can be tuned in either case. Electric dipoles 
can be polarized simply by immersing them in an elec- 
tric field, whereas magnetic dipoles may by tuned by the 
scheme described in Ref. [2^. Consequently, the ratio 
D^/b, and hence the relative importance of the dipolar 
interaction compared to the short-range interaction, can 
be changed essentially at will. This motivates us to in- 
vestigate the zero temperature equilibrium properties of 
dipolar Bose gases over a wide range of I?*/6, includ- 
ing the short-range dominated regime with D<,/6 <C 1 
and the long-range dominated regime with £)*/6 ;3> 1. 
Note that, to date, the effects of dipolar interactions have 
been observed experimentally only for atomic Cr with 
D^/b « 0.2 (taking b to be the s-wave scattering length). 

We start our discussion by considering two interact- 
ing dipoles, i.e., we set A^ = 2 in Eq. (^. After sepa- 
rating off the center of mass part oi H , we rewrite the 
Hamiltonian for the relative coordinate in spherical co- 
ordinates and solve the corresponding two-dimensional 
Schrodinger equation numerically using standard tech- 
niques. We first determine scattering and bound state 
solutions in the absence of an external confining poten- 
tial, i.e., for w = 0. Figure PT^a) shows the zero-energy 
s-wave scattering length o as a function of D^/b. We 
refer to a calculated for d^ = as the "bare" scattering 
length and to a calculated for finite cP as the "dipole- 
normalized" scattering length. For Z?* = 0, no two-body 
bound states exist and a is equal to b. The scattering 
length a decreases with increasing 13*, and diverges and 
changes sign at D^,/b sa 8.5, signaling the creation of a 
two-body bound state. At D^,/b « 19, a shows a second 
divergence corresponding to a second s-wave bound state 
being pulled in. 

Trapped two-body systems could be prepared experi- 
mentally by loading ultracold polar molecules into a very 
deep optical lattice and realizing doubly-occupied lattice 
sites. The two-body energy for different electric field 
strengths can then be measured spectroscopically. To de- 
termine the energy spectrum of two trapped dipoles we 
fix the short-range two-body length 6, i.e., b = 0.0137a/io, 
and vary D^^/auo, where auo denotes the oscillator length, 
aho = yjhjirfiuj). Dashed fines in Fig. QJb) show the 
total energy EjN per dipole as a function of D^,/aho- 
Comparison of Figs. H^a) and (b) reveals that the ener- 
getically lowest-lying state with positive energy becomes 
negative at about the same value of I?* as that for which 
the scattering length a diverges [the D* values shown 
in panels (a) and (b) of Fig. ^ extend, although scaled 
differently, over the same range]. Furthermore, the trap 
energies nearly coincide with those for a non-interacting 
two-particle gas, i.e., E/N = 1.5, 2.5, • • • hw, at D* values 
for which a = 0. The s-wave scattering length, which de- 
pends only on the ratio -D*/&, thus determines the gross 
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FIG. 1; (a) Solid lines show the s-wave scattering length a as 
a function of the dipole strength D, , both in units of b, for the 
two-body potential V, Eq. (|5J. Vertical dotted lines denote 
those D, /b values at which a diverges and a new bound state 
appears in the two-body potential, (b) Dashed lines show 
E/N for two dipoles under external spherical confinement cal- 
culated for b = 0.0137aho as a function of Dt/ttho obtained 
by solving the linear Schrodinger equation for the Hamilto- 
nian given by Eq. Q . Solid lines show the corresponding GP 
energy obtained by solving Eq. (|1J for the pseudo-potential 
VcB, Eq. @, using the dipole-normalized scattering length. 
The GP energies are plotted for each branch of the two-body 
spectrum. Note that the D* values shown in panels (a) and 
(b) extend over the same range. 



features of the energy level spectrum of two interacting 
dipoles under external spherical confinement. The details 
of the energy spectrum such as the slope of the energy 
levels near E/N = {n + l/2)fiuj, however, depend addi- 
tionally on the magnitude of D^ or b. 

For iV > 2, we solve the Schrodinger equation using the 
DMC technique with importance sampling, which deter- 
mines the ground state energy of the time-independent 
Schrodinger equation by propagating an initial "walker 
distribution" in imaginary time and projecting out the 
lowest stationary eigenstate |25j . To efficiently treat large 
systems, a stochastic realization of the short-time Green's 
function propagator is used, which introduces a statis- 
tical uncertainty of any DMC expectation value. De- 
tails of the procedure will be presented elsewhere. Sym- 
bols in Fig. 121 show our DMC energies E/N per dipole 
for b = 0.0137a/io (and cP values for which V supports 
no two-body bound states) as a function of D^/aho for 
A^ = 4, 10, 20 and 50. Statistical uncertainties are indi- 
cated by vertical error bars. For completeness, dashed 
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FIG. 2: Symbols show the energy per particle E/N calcu- 
lated by the DMC method for b = 0.0137aho as a function of 
Dt/aho for N = 4, 10, 20 and 50. Vertical error bars indi- 
cate statistical uncertainties. For completeness, a dashed line 
shows E/N, calculated using B-splines, for N = 2. Solid lines 
show E/N calculated by solving the non-hnear GP equation, 
Eq. Q, with the dipole-normalized scattering length. For 
comparison, a dotted line shows E/N for TV = 10 calculated 
by solving the non-linear GP equation, Eq. I^J i with the bare 
scattering length, i.e., the cut-off radius b. 



lines show the E/N data for TV = 2 from Fig. ^h). 
The energy E/N per dipole decreases with increasing 
D^. In particular, E/N becomes smaller than the ideal 
gas value of l.Sfio; for negative s-wave scattering lengths 
(D*/aho greater than w 0.06 in the figure). Finally, for 
fixed D^/ttho the attractive part of the dipolar interaction 
leads to a decrease of E/N with increasing N. We find 
qualitatively similar behaviors for dipolar gases confined 
in elongated cigar-shaped and pancake-shaped traps. 

Our variational many-body calculations for dipolar 
gases show that the region in configuration space where 
the metastable condensate exists is separated by an "en- 
ergy barrier" from the region where bound many-body 
states exist. This energy barrier is familiar from varia- 
tional treatments of atomic EEC's with attractive inter- 
actions |26|,|23,|28j. The existence of this barrier is crucial 
for our DMC calculations to converge to the metastable 
condensate state for sufficiently large D^/aho and not to 
the cluster-like ground state. The dipolar gas collapses at 
the D^/aho value for which the energy barrier vanishes. 
Our DMC calculations show that the condensate prior 
to collapse is only slightly elongated, which is consistent 
with our finding that the collapse is induced primarily by 
the negative value of a. 

We now assess the validity of the GP equation for 
trapped dipolar Bose gases, which can be derived by per- 
forming a functional variation of the expectation value 
of the Hamiltonian given by Eq. (Q , calculated with re- 
spect to a product wave function ^P, '^(fi,--- ,rAr) = 



Y\j=iX{^j)- For this procedure to be meaningful, the 
two-body interaction potential V, Eq. (01, has to be re- 
placed by a pseudo-potential V^s l^]'- 



Kff(r) 



47rft2a(d) 



m 



,1-3 cos^ e 



(3) 



whose zero-energy T-matrix, calculated in the first Born 
approximation, reproduces the full zero-energy T-matrix 
of the model potential V , Eq. Q. The strength of the 
contact term of Vcs is not, as might be expected naively, 
given by the cutoff radius b but by the dipole-normalized 
s-wave scattering length a{d). The GP equation for the 
single particle orbital x(r) then reads 



2m 2 m 
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X{r} = ex(r): (4) 



where e denotes the chemical potential. We solve the non- 
local Eq. Q numerically by the steepest descent method. 
At each time step, the integration over the dipole poten- 
tial is evaluated in momentum space with the aid of fast 
Fourier transforms [ig. Once the solution to the GP 
equation is found for a given TV, a and d^, the total en- 
ergy E can be obtained straightforwardly. 

Solid lines in Fig. |21 show the GP energies £'/TV per 
dipole for various TV as a function of D^,/aho- Figure [21 
shows excellent agreement between the GP and DMC 
energies (symbols) for all TV considered. To illustrate 
that this agreement depends crucially on the value of the 
dipole-normalized s-wave scattering length a in the con- 
tact part of the pseudo-potential Vcs, Eq. Q, a dotted 
line in Fig.[21shows the GP energy per dipole for TV = 10 
obtained using the cutoff radius b instead of a. Figure [21 
indicates that this simple description overestimates £'/TV 
severely when the dipole length £>* becomes comparable 
to and larger than the short-range length b. 

The dipole-dependent scattering length has important 
implications for the stability of a condensate, as shown 
in the TV-vs-dipole stability diagram in Figure [21 For 
concreteness, we have included an alternative horizontal 
axis, representing the dipole moment in Debye, assuming 
a trap of frequency ;/ = 1 kHz, and a molecular mass 20 
amu, typical for light molecules. The shaded and white 
areas in Fig.[31denote parameters for which the GP equa- 
tion does and does not possess a solution, respectively. 
The dark shaded areas represent where the condensate is 
expected to be stable for any number of molecules, even 
in free space, as given by the criterion a{d) > D^/12tt 
[29j. Apart from these regions of "absolute" stability, 
the condensate for a fixed dipole moment will ultimately 
become unstable as the number of molecules is increased. 
Indeed, for certain values of dipole where a{d) takes large, 
negative values (say, near D^/a\^o = 0.27), a condensate 
is not supported at all. 
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FIG. 3: Partial stability diagram for a dipolar condensate 
The white regions depict parameters for which the gas is pre- 
dicted to be mechanically unstable. Shaded areas are regions 
of stability, and dark shaded areas denote parameters that 
are expected to produce stable condensates even in free space. 
The top axis translates the dipole length Z), into a dipole mo- 
ment, assuming a trap of frequency v = IkHz and a molecular 
mass 20 amu. 



Alternatively, for fixed A'", the condensate stability can 
be probed as a function of dipole moment. This is likely 
a parameter more amenable to fine-tuning in the labo- 
ratory. In this case, an initially stable condensate will 
collapse after the dipole exceeds a certain value. There 
then follows a region of instability, followed by another 
region of stability as the dipole is made yet larger and the 
scattering length takes positive values. This alternating 
pattern of stable and unstable condensates continues be- 
yond the two-and-a-half cycles we have shown in Fig. 13 
This pattern is in contrast to the generally held view 
of polar condensate collapse, which would posit a single 
collapse when the dipole reaches a large critical value. 
Instead, there are many critical values, generated each 
time a new bound state is absorbed into the two-body 
potential. Because this collapse is largely s-wave dom- 
inated, the gas is much more nearly isotropic near the 
collapse point than has previously been reported. 

In summary, we have tested, for the first time, the 
validity of the GP equation for describing Bose- Einstein 
condensates interacting via strong dipolar forces. We find 
that the GP equation works quite well as compared to 
essentially exact DMC methods, as long as the depen- 
dence of the s-wave scattering length on dipole moment 
is accounted for. Doing so, we predict a rich stability 
diagram for such a system, incorporating alternating re- 
gions of stability and instability as the dipole moment is 
varied. 
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